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Abstract 


Following  the  visit  by  Brad  Johnson,  Joe  Gailani  and  Ty  Wamsley  in  February  2006  and 
after  further  discussions  the  (at  least  short-term)  line  of  approach  has  been  to  adopt 
ShoreCirc  as  the  nearshore  vehicle  of  development,  and  to  investigate  the  possibilities  to 
either  extend  SC  to  include  overland  flow,  breaching  etc,  or  to  couple  it  with  a  model 
specifically  for  the  really  nearshore,  swash  and  above. 

Ap  van  Dongeren  has  since  downloaded  ShoreCirc  and  has  compiled  the  code  on  Delft 
Hydraulics’  computers.  He  has  since  tried  to  couple  the  code  with  the  existing  wave  driver 
developed  by  Ozkan,  but  at  this  point  the  coupling  still  has  bugs.  He  has  been  in  contact 
with  Brad  Johnson  at  University  of  Delaware  and  Tuba  Ozkan  about  these  problems.  He  is 
planning  a  trip  to  Delaware  in  the  coming  months,  to  discuss  problems  and  further 
development. 

Dirk-Jan  Walstra  has  been  working  with  Jaap  van  Thiel  de  Vries  on  the  problem  of 
modelling  the  dynamics  of  dune  foot  retreat,  which  has  issues  such  as  whether  the  dune 
erosion  just  follows  the  transport  of  sand  in  the  wet  profile  or  whether  it  is  necessary  to 
include  processes  to  explicitly  get  the  sand  out  of  the  dune  by  some  physical  process,  and 
how  to  model  this.  They  have  been  developing  a  schematised  Matlab  environment  to  better 
test  the  concepts,  without  the  burden  of  a  big  model  system. 

Dano  Roelvink  and  Ad  Reniers  have  been  looking  at  the  feasibility  of  developing  a  separate 
model  for  the  really  shallow,  swash  and  dune  area  and  have  made  tangible  progress.  They 
are  developing  a  set  of  Matlab  functions  for  short  wave  propagation,  instationary  shallow 
water  equations,  sediment  transport  and  continuity  equations  that  can  be  coupled  in  various 
ways  and  are  designed  to  cope  with  extreme  conditions  such  as  encountered  during 
hurricanes.  Since  length  scales  are  short  in  terms  of  wave  lengths  and  supercritical  flow  is 
happening  a  lot,  the  numerical  implementation  is  mainly  first  order  upwind,  which  in 
combination  with  a  staggered  grid  makes  it  quite  robust;  to  keep  the  code  simple  and  make 
coupling  and  parallellization  easier,  we  have  used  explicit  schemes  with  automatic  time¬ 
stepping,  thereby  guaranteeing  stability. 

The  short  wave  propagation  contains  a  newly-developed  time-dependent  version  of  HISWA, 
which  solves  the  wave  refraction  and  allows  variation  of  wave  action  in  x,y,  time  and  over 
the  directional  space,  and  can  be  used  to  simulate  the  propagation  and  dissipation  of  wave 
groups;  the  nice  thing  compared  to  the  existing  surfbeat  model  is  that  a  separate  wave  model 
to  predict  the  mean  wave  direction  is  not  needed,  and  that  it  allows  different  wave  groups  to 
travel  in  different  directions.  We  have  done  some  principle  tests  so  far  and  are  able  to 
simulate  surfbeats  running  up  and  over  dunes.  Ad  Reniers  has  included  full  wave-current 
interaction  in  the  short  wave  propagation. 

A  basic  sediment  transport  module  has  been  included  which  solves  the  2DH  advection- 
diffusion  equation  and  produces  total  transport  vectors,  which  can  be  used  to  update  the 
bathymetry.  The  pickup  function  so  far  has  been  a  very  simple  function;  in  the  coming 
period  a  more  realistic  version,  following  Reniers  et  al  (2004),  will  be  implemented. 
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We  have  looked  at  the  possibility  of  coupling  such  a  model  with  Shorecirc,  by  first  coupling 
it  with  itself,  cutting  the  domain  in  half  and  exchanging  information  two  ways.  This  goes 
without  any  problems,  so  we  guess  that  coupling  with  Shorecirc  will  be  doable  as  well. 

The  coming  period  will  be  mainly  devoted  to  testing  of  the  various  components  against 
analytical,  laboratory  and  field  test  cases  and  to  investigate  the  coupling  with  ShoreCirc. 


1  Introduction 


This  report  is  the  first  progress  report  of  the  project  ‘Modeling  of  Hurricane  Impacts’, 
contract  no.  N62558-06-C-2006,  which  was  granted  by  the  US  Army  Corps  of  Engineers, 
Engineer  Research  and  Development  Center  (ERDC),  European  Research  Office  and 
administered  by  FISC  SIGONELLA,  NAVAL  REGIONAL  CONTRACTING  DET 
LONDON,  SHORE/FLEET  TEAM. 

The  project  is  being  carried  out  by  Prof.  Dano  Roelvink  of  UNESCO-IHE  (Principal 
Investigator),  Dr.  Ad  Reniers  of  Delft  University  of  Technology  and  Dr.  Ap  van  Dongeren 
and  Dirk-Jan  Walstra  of  WL  |  Delft  Hydraulics. 

Following  the  visit  by  Brad  Johnson,  Joe  Gailani  and  Ty  Wamsley  in  February  2006  and 
after  further  discussions  the  (at  least  short-term)  line  of  approach  has  been  to  adopt 
ShoreCirc  as  the  nearshore  vehicle  of  development,  and  to  investigate  the  possibilities  to 
either  extend  SC  to  include  overland  flow,  breaching  etc,  or  to  couple  it  with  a  model 
specifically  for  the  really  nearshore,  swash  and  above.  These  activities,  which  fall  under 
tasks  4  and  7  of  the  project  proposal,  were  indicated  as  having  top  priority  and  have 
consequently  received  most  attention  in  the  past  3  months. 

In  the  following  we  will  outline  in  some  detail  the  activities  and  results  over  this  period.  In 
Chapter  2  the  activities  related  to  ShoreCirc  are  explained.  Chapter  3  gives  an  overview  of 
the  analysis  of  dune  erosion  mechanisms  and  their  relation  to  nearshore  profile 
development.  Chapter  4  gives  a  description  of  the  new  model  with  working  title  XBeach. 
Chapter  5  outlines  the  proposed  activities  for  the  coming  3  months. 


2  ShoreCirc  adaptations 


The  Shorecirc  model  has  been  developed  by  the  research  group  of  the  late  lb  Svendsen  of 
the  University  of  Delaware.  Dr.  Ap  van  Dongeren  was  part  of  that  research  group  and  is  one 
of  the  co-writers  of  the  code.  At  Delft,  he  has  implemented  the  time-dependent  energy  and 
roller  equations  (Van  Dongeren  et  al.,  2003)  into  the  code.  These  extensions  are  not  yet  in 
the  motherversion  and  will  be  implemented  in  this  project  in  the  next  quarter.  In  fact,  the 
extension  of  the  energy  equation  as  described  above  will  be  implemented. 
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Next  the  Shorecirc  code  (for  the  midfield  of  a  few  surfzone  widths)  will  be  coupled  to  the 
Xbeach  model  (very  nearshore  and  dune)  and  to  output  of  the  SWAN  model  (far  field  model 
for  quasi-linear  wave  propagation  over  areas  of  order  100  km). 

3  Dune  erosion  mechanisms 


3.1  Introduction 

Within  the  Morphos  project  WL  |  Delft  Hydraulics  and  TU  Delft  are  working  on  the 
modeling  of  dune  erosion.  Although  many  empirical  formulations/models  exist  we  are 
focusing  on  process  based  approaches.  The  problem  with  dune  erosion  is  that  the  modeling 
of  the  retreat  of  the  upper  beach  and  dune  front  requires  additional  formulations  which 
cannot  be  based  on  the  traditional  transport  modeling  techniques.  Various  approaches  have 
been  developed  over  the  years  which  are  primarily  based  on  an  extrapolation  of  the  offshore 
directed  sediment  transport  in  some  “last  wet”  point  (see  e.g.  Steetzel,  1993).  The  most 
advanced  process  based  model  that  works  with  such  an  approach  is  the  DurosTa  model 
(Steetzel,  1993).  This  model  has  undergone  extensive  validation  and  has  proven  to  yield 
reliable  and  robust  predictions  for  a  range  of  conditions.  However,  this  model  contains 
several  novelties  of  which  their  importance  to  an  accurate  dune  erosion  retreat  is  unkown. 
Issues  that  require  further  investigation  are  the  extrapolation  method  and  the  transport 
formulations.  Here  we  will  be  discussing  our  findings  on  the  extrapolation  methods. 


3.2  Extrapolation  methods 

As  a  first  step  we  investigated  the  effect  of  various  extrapolation  methods.  All  extrapolation 
methods  are  based  on  a  preset  transport,  Sbase,  which  is  constant  in  time.  In  total  three 
methods  were  investigated.  The  first  method  is  developed  in  this  project  and  involves  a 
direct  scaling  of  the  transport  using  the  local  water  depth. 

=  -j-1  (o.i) 

y  base  ) 

The  second  method  uses  a  vertical  extrapolation  which  results  in  a  dune  erosion  retreat 
where  the  shape  of  the  dune  is  maintained. 

r  /  \  A 


5,  (x)  =  S, 


base 


1 


'base 


Z  —  Z 
top  base  J 


(0.2) 


in  which  z  is  the  local  profile,  subscripts  base :  point  for  last  wet  point,  top\  indicates  the 
position  of  the  dune  top. 


The  third  method  uses  the  DurosTa  formulation  (see  Steetzel,  1993)  for  details. 


Some  prilimenary  results  for  these  extrapolation  methods  are  shown  in  Figure  1  for  a 
schematic  profile  in  which  the  base  transport  was  kept  constant  at  1  m3/m3.  In  this  figure  the 
profile  development  (top  row),  sedimentation-erosion  (middle  row)  and  the  cross-shore 
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sediment  transport  (bottom  row)  are  shown  in  time.  The  bold  lines  indicate  the  initial 
distributions. 


What  becomes  apparent  from  these  results  is  that  method  1  (left  column)  yields  very  similar 
predictions  as  the  DurosTa  method  (right  column).  The  steepening  of  the  dune  face  is 
represented  in  both  methods,  whereas  method  2  (middle  column)  results  in  a  dune  retreat 
with  a  constant  shape. 
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Figure  1  Dune  erosion  according  to  three  different  extrapolation  methods.  Top  row: 
development  of  dune  profde,  second  row:  bottom  changes  per  time  step,  third  row:  cross¬ 
shore  sediment  transport  across  dune.  Left  column:  Method  1,  middle  column:  Method  2, 
right  column:  Method  3. 


3 .3  Anticipated  work 

At  the  moment  we  are  extending  the  work  presented  above  by  looking  into  the  required 
numerical  schemes  and  stability.  Furthermore  a  detailed  comparison  between  method  1  and 
DurosTa  is  required  to  assess  the  necessity  of  including  many  of  the  features  that  are 
available  in  the  DurosTa  model.  In  addition  we  will  investigate  various  transport 
formulations  and  the  effects  of  surfbeat  on  these  types  of  extrapolation  techniques. 
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4  XBeach  model  development 


4.1  Objectives 

The  main  objective  of  the  XBeach  model  is  to  provide  a  robust  and  flexible  environment  in 
which  to  test  morphological  modelling  concepts  for  the  case  of  dune  erosion,  overwashing 
and  breaching.  The  top  priority  is  to  provide  numerical  stability;  first  order  accuracy  is 
accepted  since  there  is  a  need  for  small  space  steps  and  time  steps  anyway,  to  represent  the 
strong  gradients  in  space  and  time  in  the  nearshore  and  swash  zone.  Because  of  the  many 
shock-like  features  in  both  hydrodynamics  and  moiphodynamics  we  choose  upwind 
schematizations  as  a  means  to  avoid  numerical  oscillations  which  can  be  deadly  in  shallow 
areas. 

The  modelling  environment  should  be  flexible  and  the  code  easy  to  comprehend  and 
concise;  therefore  we  have  adapted  the  Matlab  environment  as  development  environment; 
conversion  to  Fortran  90/95  at  a  later  stage  will  be  straightforward. 

The  code  should  have  the  following  functionality: 

•  Depth-averaged  shallow  water  equations  including  time-varying  wave  forcing 
terms;  combination  of  sub-  and  supercritical  flows; 

•  Time- varying  wave  action  balance  including  refraction,  shoaling,  current  refraction 
and  wave  breaking; 

•  Wave  amplitude  effects  on  wave  celerity; 

•  Depth-averaged  advection-diffusion  equation  to  solve  suspended  transport; 

•  Bed  updating  algorithm  including  possibility  of  avalanching; 

•  Possibility  of  online  coupling  with  ShoreCirc; 

•  Possibility  to  extend  to  parallel  multi-domain  version; 

4.2  General  setup 

The  program  XBeach  consists  of  a  main  Matlab  script,  xbeach.m,  and  a  number  of functions 
that  operate  on  two  structures'. 

•  par  -  this  contains  general  input  parameters 

•  s  -  this  contains  all  the  arrays  for  a  given  computational  domain 

For  a  single-domain  run,  one  structure  s  is  passed  between  flow,  wave,  sediment  and  bed 
update  solvers,  which  extract  the  arrays  they  need  from  the  structure  elements  to  local 
variables,  do  their  thing  and  pass  the  results  back  to  the  relevant  structure  elements.  This 
makes  the  overall  program  clear,  prevents  long  parameter  lists  and  makes  it  easy  to  add 
input  variables  or  arrays  where  needed. 

For  multi-domain  runs,  one  can  define  multiple  instantiations  of  the  structure  s  which  are 
passed  to  the  same  functions;  an  additional  function  is  needed  to  pass  the  boundary 
information  between  the  domains  back  and  forth.  We  have  carried  out  a  simple  test  of  this 
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principle,  without  actually  implementing  a  multi-processor  version,  which  confirms  that  the 
data  structure  can  handle  this  case. 

In  the  Table  1  we  will  outline  the  various  functions  and  their  purposes.  The  main  script 
xbeach.m  is  reproduced  in  Table  2  below. 


Function  call 

Purpose 

[par]  =  waveinput 

Creates  elements  of  structure  par  containing  wave  input 
parameters 

[par]  =  flowinput(par) 

Adds  elements  of  structure  par  containing  flow  input 
parameters 

[par]  =  sedinput(par) 

Adds  elements  of  structure  par  containing  sediment  input 
parameters 

[s]  =  gridbathy; 

Creates  grid  and  bathymetry  and  stores  them  in  structure  s 

[s]  =  wave_dist  (s,par); 

Creates  initial  directional  spectrum  at  sea  boundary 

[s,par]=wave_init  (s,par); 

Initialises  arrays  (elements  of  s)  for  wave  computations 

[s]  =  flow_init  (s,par); 

Initialises  arrays  (elements  of  s)  for  flow  computations 

[s]  =  sed_init  (s); 

Initialises  arrays  (elements  of  s)  for  sediment  computations 

[s]  =  wave_bc  (s,par,it); 

Wave  boundary  conditions  update,  each  timestep 

[s]  =  flow_bc  (s,par,it); 

Flow  boundary  conditions  update,  each  timestep 

[s]=wave_timestep  (s,par); 

Carries  out  one  wave  timestep 

[s]=flow_timestep  (s,par); 

Carries  out  one  flow  timestep 

[s]=transus(s,par); 

Carries  out  one  suspended  transport  timestep 

[s]=bed_update(s,par); 

Carries  out  one  bed  level  update  timestep 

output(it,s,par) 

Performs  output 

Table  1  Overview  of  Matlab  functions  XBeach 


Of  these  functions,  the  input  functions  and  grid  bathy  are  case-dependent  at  present  and 
actually  contain  the  input  values  and  grid  and  bathymetry  definitions.  This  can  however 
easily  be  replaced  by  functions  with  the  same  output,  which  read  data  from  input  files  or 
input  them  through  a  screen  dialog.  We  aim  to  leave  the  initialisation,  boundary  conditions 
and  computational  functions  case-independent.  The  output  function  can  be  adapted  to  fit 
specific  needs. 
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clear  all 

"o 

%  General  input  per  module 
[par]  =  wave_input; 

[par]  =  f low_input (par ) ; 

[par]  =  sed_input  (par ) ; 

%  Grid  and  bathymetry 
[s]  =  grid_bathy; 

%  Directional  distribution  wave  energy 
[s]  =  wave_dist  (s,par); 

%  Initialisations 

[s,par]  =  wave_init  (s,par); 

[s]  =  flow_init  (s,par); 

[s]  =  sed_init  (s); 

nt=par . nt ; 

for  it=l:nt; 

"o 

%  Wave  boundary  conditions 
[s]  =  wave_bc  (s, par, it); 

%  Flow  boundary  conditions 
[s]  =  flow_bc  (s, par, it); 

%  Wave  timestep 

[s]  =  wave_timestep  (s,par); 

%  Flow  timestep 

[s]  =  f low_timestep  (s,par); 

%  Suspended  transport 
[ s ] =transus ( s , par ) ; 

%  Bed  level  update 
[s] =bed_update (s,par) ; 

%  Output 
output ( it , s , par ) 

end 

Table  2  Main  routine  xbeach.m 


4.3  Matlab  code  and  test  case 

In  the  attached  xbeachv001.zip  the  Matlab  code  is  provided,  for  a  test  case  of  a  plane  sloping 
beach  with  a  Gaussian  hump  superimposed,  and  longcrested  wave  groups  are  entering  at  a 
20  degree  angle  to  the  shore  normal.  The  grid,  beach  slope  and  hump  dimensions  can  be 
easily  modified  in  the  function  grid  bathy.m  and  incident  wave  conditions  that  drive  this 
case  are  defined  in  wave  input.m. 

In  Figure  2  below  a  snapshot  of  the  water  level  and  velocity  pattern  is  shown.  The  model 
can  clearly  handle  wetting  and  drying  and  (breaking)  long  waves. 

The  details  of  the  most  important  elements,  the  wave  action  equation  solver  and  the  shallow 
water  equations  solver,  are  given  in  the  next  sections. 


Netherlands  Centre  for  Coastal  Research 


4-7 


Modeling  of  hurricane  impacts 
Progress  Report  1  -  March-May,  2006 


1000 

900 

800 

700 

600 

500 

400 

300 

200 

100 

0 

Figure  2  Snapshot  of  water  level  and  velocity,  beach-with-hump  testcase. 


4.4  Wave  action  equation  solver 

The  wave  forcing  in  the  shallow  water  momentum  equation  is  obtained  from  a  time 
dependent  version  of  the  wave  action  balance  equation.  Similar  to  Delft  University’s 
HISWA  model,  the  directional  distribution  of  the  action  density  is  taken  into  account 
whereas  the  frequency  spectrum  is  represented  by  a  single  mean  frequency.  The  wave  action 
balance  is  then  given  by: 

8A  dcA  8c  A  8c  „  A  D 
—  +  — —  +  — —  +  — —  = - 

dt  dx  8y  86  a 

with  the  wave  action: 


A(x,  y,  6) 


Sw(x,y,0) 

o-(x,y) 


where  Sw  represents  the  wave  energy  in  each  directional  bin  and  <J  the  intrinsic  wave 
frequency.  The  wave  action  propagation  speeds  in  x-  and  y-direction  are  given  by: 

cx(x, y, 6)  =  cg (x,  y)-cos(6)  +  u(x, y) 
cy  (x,  y,  6)  =  cg  (x,  v>sin(<9)  +  v(x,  y) 
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This  concludes  the  advection  of  wave  action.  The  wave  energy  dissipation  due  to  wave 
breaking  is  modelled  according  to  Baldock  et  al.  [1998]: 

D  =  ^aQlPgf„(Hl+Hl) 

with  a  =  (9(1)  and  fm  representing  the  mean  intrinsic  frequency.  The  fraction  of  breaking 
waves  is  given  by: 
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Qb  =  exP 


1 

j_r 2  I 

V  1  rms  J 

where  the  breaking  wave  height  is: 
ykh 


a88tanh 


0.88 


and  ^  is  a  calibration  parameter.  The  root  mean  square  wave  height  is  obtained  from: 


TT 

rms  ^ 


8  fsw(x,y,0)d0 


Pg 


Next  the  total  wave  dissipation,!),  is  distributed  proportionally  over  the  wave  directions: 


D(x,y,&) 


Sw(x,y,0)- 

Ew(x,y) 


This  closes  the  set  of  equations  for  the  wave  action  balance.  Given  the  spatial  distribution  of 
the  wave  action  and  therefore  wave  energy  the  wave  forcing  can  be  calculated  utilizing  the 


radiation  stress  tensor: 

F  =- 

fas„  ,  ^ 

X 

v  dx  dy  , 

F  =  - 

(M*  ,  ^ 

y 

v  dx  dy  j 

And: 

( 

(  r 

*  ^  rr  /  O 

x  1^1 

S„  = 

— (l  +  cos 

0)-- 

XX  J 

v 

>  2 
z  J 

^ xy  =  E yx  —  j"  sin  0  COS  0 


SwdO 

\ 

J 


\d0 


\{—  (l  +  sin2  6 

)--) 

>  2 
zy 

sje 


Similar  to  the  solution  of  the  shallow  water  equations  we  use  an  up-wind  schematisation  to 
solve  the  wave  action  balance.  The  wave  action  is  given  at  the  same  points  at  the  water 
level.  The  advection  of  wave  action  is  then  discretized  as  follows: 
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dcnxA" 

dx 

dc”A” 

dx 

dc"A" 

y 

dc"  A" 

y 

dy 


(i,j,k)  = 


(ij,k)  = 


(ij,k)  = 


(ij,k)  = 


rn  An  —rn  An 

L'x,i,j,k^Li,j,k  L'x,i-lJ,kjrLi-lJ,k 

x,rx^ 

rn  A"  —r"  An 

{-'x,i+\,j,k  M,j,k  i,jX 

xi+ij~xio 

rn  A "  -r"  A" 

{-'y,i,j,k  i,j,k  ^y,i,j-\,k  i,]-\,k 

J'g-J'y-i 

rn  An  —rn  An 

^y,i,j+\,k  Uj+l.k  L'y,i,j,ki,j,k 


y-m-y* 


n  j  n  n  a  n  n  a  n 

dc0A  .  ce,i,j,kA,j,k  ~  ce,i,j.k-\Ai  jk-\ 

-0,J,k)  =  — 


do 


~  n  jn  n  a  n  n  ah 

dc0A  .  Cg,i,j,k+l^iJ,k+l~C0,i,j,kA,j,k 

(i,J,k)  =  — 


do 


#y,k+i  -^ij,k 


’  Cx,i,j,k 


’  Cx,i,j,k 


’  Cy,i,j,k 


’  Cy,iJ,k 


’ C0,i,j,k 


’  C0,i,j,k 


>0 

<0 

>0 

<0 

>0 

<0 


Similar  for  the  wave  action  balance: 


A”+l  _  An 

Ai,],k  Ai,j,k 

At 


dc^i 

dx  ij,k 


d<A 

' 


i.j,k 


dcngAn  D 

dd  ij,k  (J  i,j,k 


which  yields  the  wave  energy  at  the  new  time  level. 


4.5  Shallow  water  equations  solver 


Shallow  water  equations,  neglecting  Coriolis  and  horizontal  diffusion  terms,  and  (grey 
terms),  for  the  moment,  wind  shear  stress: 


du  du  du  r  r, 

—  +  u  —  +  v —  =  — 
dt  dx  dy  ph  ph 

dv  dv  dv  t  zb 

- h  U - b  V - =  +— - — 

dt  dx  dy  ph  ph 

dr)  dhu  dhv 

— --f - -i - = 

dt  dx  dy 

Here,  h  is  the  water  depth,  u,  v  are  velocities  inx  and  and  y  direction,  zbx,  rhv  are  the  bed 


dti  F 
g  —  +  — 
dx  ph 

(2.1) 

-e’F  +  F 

dy  ph 

(2.2) 

) 

(2.3) 

shear  stresses,  g  is  the  acceleration  of  gravity,  r/  is  the  water  level  and  F  ,  Fr  are  the 
wave-induced  stresses. 


We  apply  an  upwind  schematisation,  since  the  horizontal  scale  of  the  problem  is  limited  and 
such  a  scheme  deals  with  shocks  in  a  natural  way. 
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We  apply  a  staggered  grid,  where  bed  levels  and  water  levels  are  defined  in  the  centre  of 
cells,  and  velocity  components  at  the  cell  interfaces. 

If  nx,ny  are  the  number  of  cells  in  both  directions,  the  water  level  points  are  numbered  from 
1  to  nx+1  and  from  1  to  ny+1. 

The  water  level  gradients  are  computed  at  the  cell  interfaces  and  are  given  by: 


dr/ 

dx 

d?i 

¥ 


II 

+ 

1 

(2.4) 

II 

+ 

1 

CL-. 

(2.5) 

For  computing  the  shear  stresses  at  the  cell  interfaces  we  need  the  velocity  magnitudes  at 
these  interfaces.  These  are  composed  by  combining  the  normal  velocity  component  at  the 
interface  and  the  average  of  the  4  adjacent  tangential  components: 


]_ 

4 

1 


vu,ij  =  T  (v;  j-i  +  vi,j  +  VMJ- 1  +  vMj) 


(2.6) 


Uv,iJ  = +  U,J  +  WMJ+1  +  UiJ+l) 


The  water  depth  in  each  cell  is  computed  as: 


hu  =  tlu-zb 


(2.7) 


The  depth  at  the  interfaces  is  taken  as  the  upwind  depth  in  case  the  velocity  is  greater  than  a 
minimum  velocity,  or  the  average  depth  between  the  adjacent  cells  in  case  the  velocity  is 
less  than  this  minimum  velocity: 


h 

=  h 

u.  . 

>  U- 

U,l,J 

hj 

hj 

min 

K,ij 

=  K,i 

'  UUJ 

<  —u 

min 

(2.8) 

K,ij 

\ui 

\<U  • 

1 1  min 

h  . . 

=  h 

v.  . 

>  v  ■ 

V,l,J 

*>J  5 

hj 

min 

Ktj 

=  h 

niJ+\  5 

ViJ 

<  -Vmin 

(2.9) 

K, 

=\(hj+hj+ 1)’ 

\VU 

<  Vmi„ 

The  advection  terms  in  x-direction  are  approximated  as  follows: 

~  n  n  n 

HU  U  —  U  , 

*"  ,J  “ l’J  ,ulj> 0 


*i-ij 

u  —  =  u.  .  — - - - 

dxij  'J  x" :  —x" 


>,j 


i-hj 


du" 


u  ,,  —  U  ■ 

n  i+\,j  i,j 

u  —  =  U-  .  — - - — 

n  lJ  n  n 

OXii  xt„  ,-x,t 


,U*<  0 


M,j 


(2.10) 


hj 


du’ 

V 


=  v , 


i,j 


u,j+l  -  u .y_! 

in  n 

y  i,j+ 1  y  i,j- 1 


The  advection  terms  in  y-direction  are  approximated  as  follows: 


(2.11) 
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,5v"  =v- 

dyu  iJK 


UJ 

avH 

fyt 


'j 


■y, 


i,j- 1 


=  V,. 


vu*~vu 


l,J 


T,'. 


7+1 


•A 


’vij  >  o 


(2.12) 


3v  „  V,.+1  .  -  Vi+U 

it  —  =  if  „ , , - - - - 

a*  J 


n»7 


(2.13) 


The  momentum  equation  is  discretized  as  follows: 

du"  du"  gui.i\ftf. 


n+ 1  « 


At 


-it 


«  2 


'/i  i.j  Vi, 


n+1  « 

v  —  v. 

l,J  l,J 


At 


=  —v 


Ox  tj 

> 

i 

h".  C2 

U,l,J 

1 

Grq 

X 

+ 

i - — 

pKu 

dv” 

dv" 

-u  —  - 

n  1  n  2  .  n  2 
8ViJ\JUv,iJ  +Vi,j 

7y+i 

g 

F 

1  y,hj 

dy 


dxi 


Tij+i  -  Tij  PK 


F  . 

(2.14) 


(2.15) 


From  this,  the  velocities  at  the  new  time  step  level  are  computed.  The  water  level  is  then 
updated  by: 

„^n+ 1  „_/2  _  n+l i  n  .  n+ 1  in  ...n+li.  „.»+ 1  r  k 

77.  .  —  77  .  U  n  ■  —  U  ,  rl  x  ■  V.  .  72  .  —  V-  .  ,/Z  .  , 

*1>J  *1>J  _  _  FJ  lJ _ /-!./  /-!,/ _ / ,j_  i j _ /,/-!  /,/-!  ^  ^ 


At 


^"Vij 


Tv 


v,ij-l 


4.6  Conclusions  and  next  steps 

The  Xbeach  model  is  an  as-simple-as-possible  testing  ground  for  algorithms  related  to  wave 
propagation,  shallow  water  flow  including  surfbeats,  wetting  and  drying,  sediment  transport 
and  morphological  updating  on  shallow  beaches  and  over  dunes.  With  some  Matlab 
experience  it  is  easy  to  define  grid,  bathymetry  and  input  conditions.  Though  the  numerical 
scheme  is  rather  simple  and  could  undoubtedly  be  improved  on  specific  points,  it  is  quite 
robust  and  seems  adequate  for  the  moment. 

In  the  coming  months  we  will  test  the  various  parts  of  Xbeach  against  analytical,  lab  and 
field  cases,  which  will  be  determined  in  consultation  with  the  ERDC  staff  involved. 
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